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Abstract —We propose a Gauss-Legendre quadrature based 
sampling on the rotation group for the representation of a band- 
limited signal such that the Fourier transform (FT) of a signal can 
be exactly computed from its samples. Our figure of merit is the 
sampling efficiency, which is defined as a ratio of the degrees of 
freedom required to represent a band-limited signal in harmonic 
domain to the number of samples required to accurately compute 
the FT. The proposed sampling scheme is asymptotically as 
efficient as the most efficient scheme developed very recently. 
For the computation of FT and inverse FT, we also develop 
fast algorithms of complexity similar to the complexity attained 
by the fast algorithms for the existing sampling schemes. The 
developed algorithms are stable, accurate and do not have any 
pre-computation requirements. We also analyse the computation 
time and numerical accuracy of the proposed algorithms and 
show, through numerical experiments, that the proposed Fourier 
transforms are accurate with errors on the order of numerical 
precision. 

Index Terms —rotation group, SO(3), sampling, band-limited 
signals, Fourier transform. 

I. Introduction 

The processing and analysis of signals defined on the 
rotation group, denoted by SO (3), in the harmonic domain 
is crucial in applications found in various fields of science 
and engineering (e.g., [l]-[7]). The transformation of a signal 
on SO(3) to the harmonic domain is enabled by the Fourier 
transform (FT) on the rotation group, which is, sometimes, 
also termed as Wigner-D transform in the literature [8], 
[9]. Clearly, the ability to compute the FT on the rotation 
group from a finite number of samples of the signal is of 
significant importance. In this work, we focus on the reduction 
in the number of samples required for the representation of a 
band-limited signal defined on SO(3) such that the FT can 
be exactly and efficiently (in terms of computational effort) 
determined from the samples of the signal. 

The design of sampling schemes on the rotation group and 
the development of computationally efficient FT’s on the rota¬ 
tion group have been actively investigated in the literature [8]- 
[15]. Among the developments in the literature, we focus 
on the sampling schemes that form regular or equiangular 
grid of samples on the rotation group and support accurate 
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computation of the FT of the signal band-limited to degree 
L (formally defined in Section II-C) [8]-[10]. 

An exact method for the computation of the FT, based on 
a sampling theorem, has been proposed for an equiangular 
sampling scheme composed of 8L^ samples with 2L samples 
along each of the three Euler angles of SO(3), [8]. Using 
the separation of variables, they also developed fast algorithm 
of complexity 0{L'^) for the computation of FT. They also 
point out that the algorithm of complexity 0{L^ log^ L) can 
also be developed as a variant of the algorithm in [16]. 
However, the algorithm of complexity 0{L^ log^ L) requires 
pre-computation of 0{L^) and is therefore only practical for 
smaller band-limits (L < 256) due to the storage require¬ 
ments [8]. 

By expanding the Wigner-d functions into Chebyshev poly¬ 
nomials and using the fast polynomial transform [17], an 
algorithm for the computation of the inverse FT of complexity 
0{L^) has been developed. The inverse FT algorithm with 
reduced complexity 0{L^ log^ L) has also been developed, 
which is shown to suffer from stability issues and therefore 
requires additional stabilization steps [18]-[20], making the 
exact complexity of the algorithm difficult to determine. For 
non-equispaced samples, the FT as an adjoint of the IFT has 
also been developed [9]. However, it is applied, using the 
Clenshaw-Curtis quadrature rule, to the equiangular sampling 
scheme that consists of number of samples. 

The sampling efficiency, defined as a ratio of number of 
Fourier coefficients required to represent a band-limited signal 
in the Fourier domain to the number of samples required to 
accurately compute FT, is the fundamental property of any 
sampling scheme. In addition to the sampling efficiency, the 
computational complexity and accuracy of the FT associated 
with the sampling scheme are of significant importance. The 
existing sampling configurations that support exact (or suffi¬ 
ciently accurate) FT attain the sampling efficiency of one third 
as we show later in the paper. Furthermore, the algorithms 
associated with the existing sampling schemes, which are 
stable and accurately compute the FT and do not have any pre- 
computation requirements, have the computational complexity 
0{L% 

We summarise the contributions of this work as follows. 

• We propose a Gauss-Legendre (GL) quadrature based 
sampling scheme on the rotation group for the exact 
computation of the FT. The proposed sampling scheme 
is asymptotically as efficient as the equally most efficient 
sampling scheme developed very recently.' 

*As an alternative to the proposed Gauss-Legendre quadrature based 
scheme, an equiangular sampling scheme, based on a sampling theorem, on 
the rotation group has been recently proposed [21] with the same asymptotic 
sampling efficiency that is attained by the proposed sampling scheme. 
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• For the proposed sampling scheme, we develop fast 
algorithms of complexity 0{L'^) for the computation of 
the FT and IFT. The developed algorithms are stable, 
accurate and do not require any pre-computation. 

• We show, through numerical experiments, that FT can be 
computed to near machine precision accuracy. 

We organize the rest of the paper as follows. The math¬ 
ematical preliminaries for the signals on the rotation group 
are summarised in Section II. In Section III, we formulate 
the FT on the rotation group, present the proposed sampling 
scheme and develop fast algorithm for the computation of 
Fourier transform. The computation time and the numerical 
accuracy are analysed in Section IV. Finally, the conclusions 
are made in Section V. 

II. Preliminaries 

A. Signals on the Rotation Group 

We consider complex-valued functions dehned 

on a rotation group SO(3) parameterized by Euler angles 
where (p G [0, 27r), '& G [0, tt] and uj G [0, 27r). 
The set of square integrable complex-valued functions dehned 
on SO(3) forms a Hilbert space L^(SO(3)) with the inner 
product for two functions / and h dehned on SO(3): 

{f,h)= I f {(fijilijOj) h{ip, djUj) dip sin (I'd doj, (1) 
Jso(3) 

where (•) denotes the complex conjugate and the integral is a 
triple integral over all rotations (p,i9,uj) [8], [9], [12], [22]. 
The inner product in (1) induces a norm ||/|| = (/, The 
functions with hnite induced norm are referred to as signals 
on the rotation group or signals for short. 

B. Wigner-D Functions 

The Wigner-79 function, denoted by i?, w), is de¬ 

hned for degree £ and orders \m\, |n| < £ as [8], [9], [22] 

, ( 2 ) 

where df^ „{'d) is the Wigner-d function [22]. We use the 
following decomposition of Wigner-d function „ (d) in 
terms of complex exponentials [10] 


Remark 1 (On the polynomial degree of Wigner-d functions): 
Since cos{m'd) = Tm'(cosd) and sin(m'd) = 
sind t7m'-i(cosd), where T^'icosid) and {/^'(cosd) 
represent the Chebyshev polynomials of degree m' of the 
hrst and second kind, respectively, it can be inferred from (4) 
that d^ „(d) is a polynomial of degree £ in cosd. 


C. Fourier Transform (FT) on the Rotation Group 

The Wigner-d? functions for integer degree £ > 0 and 
integer orders |to|, |n| < £ form a complete set of orthogonal 
functions for the space L^(SO(3)). The orthogonality relation 
for Wigner-d? functions is given by 

where dr.r' is the Kronecker delta. By completeness of the 
Wigner-d? functions, any function / G L^(SO(3)) can be 
expressed as 

OO £ t 

f{p,d,Uj)=J2 (6) 

£—0 m——£ n——£ 


where 

{fY = (fD^ ) (1) 

\J )m,n 37f2 m,n/i v' / 

denotes the Fourier coefficient of degree £, orders m and n 
and forms the Fourier (or spectral) domain representation for 
signals dehned on SO(3). 

The signal / is said to be band-limited at degree L if 
(/)m n = 0’ ^ Th® of signals band-limited at degree 

L form a subspace of L^(SO(3)), which is denoted by T-Ll- 
Let 


d^ _ ^ ^ ^ I _ — 1)(2L -I- 1) 

£—f)m——£n——£ ^ 


( 8 ) 


denote the dimension of the subspace Hl- Then dz, also 
quantihes the degrees of freedom required in the Fourier 
domain to represent a signal band-limited at L. We refer to the 
computation of Fourier coefficients {f)m n from the signal /, 
given in (7), as the Fourier transform (FT) of /. The synthesis 
of a signal / from its Fourier coefficients (/)^ given in (6), 
is referred as inverse Fourier transform (IFT). 


d^ (d) = z"”"* 


'i' 'd 


(3) 


m' ——£ 


where A^ „ = d^^YFl^Y Using the symmetry relations 
exhibited by Wigner-d functions [22], the Wigner-d function 
in (3) can be expressed as 


dL,n{^) 

where 


i 

(4) 


m'—O 


{ 1, m' = 0, 

2cos(m'd), m' Y 0, m — n even, 
2i sinfm'd), m' f 0, m — n odd. 


D. Sampling on the Rotation Group 


In order to accurately determine the FT of a band-limited 
signal on the rotation group from its samples, an equiangular 
sampling scheme, denoted by .hz, has been proposed in [8], 
which is composed of 8L^ samples on a regular grid formed 
by 2L number of samples along each Euler angle in {p, d, w). 

We dehne the asymptotic sampling efficiency, denoted by 
££(•), of any sampling scheme as a ratio of the dimension dz 
of the subspace Hl formed by the band-limited signals to the 
number of samples, denoted by N, required to compute ET of 
a band-limited signal / GdfzasL—J^oo. Since the sampling 
scheme Kl requires A(.Sz) = number of samples, it 
attains an asymptotic sampling efficiency 


E{Rl) 


dL 


lim , ^^ , 
L—>oo NfM-Y 


1 

6’ 


(9) 
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We note that the sampling scheme that also requires 8L^ 
number of samples has also been proposed in [9] for the com¬ 
putation of FT. By applying the periodic extension approach 
along "d, originally proposed in [23], the sampling scheme 
with asymptotic efficiency has been proposed just 

recently in [21]. Furthermore, the algorithms developed for 
the computation of FT and IFT for these existing schemes 
have the computational complexity O(L^). 

III. Proposed Sampling and Fourier Transeorm 

We propose a sampling scheme with an efficiency of 
2E{^l) for the representation of band-limited signals defined 
on the rotation group. We formulate the FT and IFT for 
the proposed scheme and develop fast algorithms for their 
computation. 


and L samples along i?, denoted by for integer v € 
{0,1,..., L — 1}, as the roots of the Wigner-d function 

do,oi'^) = PLicos'd) = 0, (15) 

where Pl denotes the Legendre polynomial of degree L [22]. 
Note that the sample points along as roots of the Wigner- 
d function (or Legendre polynomial) can be determined pre¬ 
cisely [24]. We refer to the proposed scheme as the Gauss- 
Legendre (GL) sampling scheme. 

C. Fourier Transform — Computation 

The FT of band-limited signal f G Hl sampled over ©^ 
can be determined by first computing Am^ni'dv) for each dy, 
following (13), as 


A. Flarmonic Formulation 


We first develop the formulation for FT and IFT. For a band- 
limited signal / G Hl, we rewrite the formulation of IFT in 
(6), by changing the order of summation and using (2), as 

L-l L-1 

f{p,d,u:)= ^ ^ (10) 

m—— (L—1) n——(L —1) 


where 

L-l 

= ( 11 ) 

1=0 

Since and are only defined for orders 

|m|, |n| < i we define (/)^^„ = = 0 for 

|m|, |n| > £ to facilitate the interchange of summations in the 
formulation of FT or IFT. 

Using the formulation of IFT in (10) and the orthogonality 
of complex exponentials, we write the FT, (7), as 
2/4-1 

= AmA'<^)dL,nW Sind dd, ( 12 ) 

where Am,n(d) can be recovered from f{ip, d, uj), again using 
(10) and noting the orthogonality of complex exponentials, as 


(*27r />27r 


A. 




f{ip,d,uj)ed^‘^+^^Uujdp. 


(13) 


B. Gauss-Legendre Sampling 

The computation of the FT requires the evaluation of 
integrals over p and w in (13) and an integral over d in 
(12). Since max(|TO|) < L and max(|n|) < L, the integrals 
involved in the computation of Ajn,n(,'d) can be replaced with 
summations if at least 2L — 1 samples of the signal are taken 
along each p and oj. Once Am,n{d) is computed, the Gauss- 
Legendre quadrature rule can be used to discretise the integral 
over d by taking only L samples along d, instead of 2L 
samples considered in [8], [9]. 

With these considerations, we propose a sampling scheme, 
denoted by ©i, that is composed of 2L—1 equiangular sample 
points along each p and uj, given by 
QiTtu 

Tu=uJu= M e {0,1,... ,2L - 2}. (14) 


Am,ni'dv) - (2L- 1)2 

2L-12L-1 

E E (16) 

u—0 w—0 

where the equivalence between (13) and (16) can be readily 
verified by substituting (10) in (13) or (16) and employing the 
orthogonality of continuous or discrete complex exponentials. 

Once Am,n{dv) in (16) is computed, the Fourier coefficient 
in {f)mn I’® determined by discretising the integral in 
(12) using the GL quadrature rule’. Since we have taken L 
sample points along d as roots of (15), we need to show that 
the integrand in (12) is a polynomial in cos 6 of degree at most 
2L — 1, so that the GL quadrature can be applied to discretise 
the integral. Following Remark 1 and the formulation for 
Am,n{'d), we note that Am,n('d), given in (11), is a polynomial 
of degree L — 1. Similarly, the integrand Am,nid)d^ ,^(d) 
in (12) is a polynomial of degree L — £ — 1, which reaches 
a maximum of 2L — 2 for £ = L — 1. Therefore, the GL 
quadrature leads to 

2/4-1 

{fYm,n = E Am,nidy)di,^„(dy)q{dy), (17) 

u=0 


where q{dv) denotes the GL quadrature weight assigned to 
each sample point dy and is given by 


q{dy) = 2(^ 


sin dy ^2 


(18) 


Remark 2 (On the efficiency of GL sampling scheme): 

The proposed GL scheme requires N{@l) = L{2L — 1)^ 
number of samples in total and thus attains an asymptotic 
sampling efficiency 

Thus, the proposed GL sampling scheme is as efficient 
(asymptotically) as the most efficient scheme developed just 
recently [21], and is more efficient, by a factor of two, than 
the schemes presented in [8] and [9]. 


^GL quadrature requires L sampling points to compute the integral exactly 
of the polynomial of degree 2L — 1. 
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D. Fourier Transform — Fast Algorithm 

Here we develop fast algorithm for the computation of FT, 
given in (17), of a band-limited signal f G FLl discretised 
over the proposed GL sampling scheme 6^. Using (3), the 
FT in (17) can be expressed (after rearranging the terms) as 


(/)L,n=*" 


. 21+1 


E /\t- At n 

‘-^m' ,n^rn,n,r. 


( 20 ) 


where Cm,n,m' = Y.t=o For the com¬ 

putation of the FT, Am^ni'^v) is first be computed using 
(16) for all |m|, |n| < L with the computational complexity 
0(L^ log L) for each and with the complexity 0{L^ log L) 
and for all -dy. Once Am,ni'&v) is computed, Cm,n,m' is then 
computed for all \m\, |n|, \m'\ < L in 0{L^), which is then 
used to compute {f)m n using (20), again with the complexity 
O(L^). Thus the overall complexity of FT is O(T^). For 
the proposed sampling scheme, we note that a fast algorithm 
of complexity 0{L^) can also be developed by using the 
separation of variables approach as adopted in [8]. 


E. Inverse Fourier Transform (IFT) — Fast Algorithm 
For the computation of IFT, we rewrite (11), using (3), as 

L-l 

AmA^y) = ^ ( 21 ) 

m' = — (L — l) 

where = J2e=o {f)L,u^L' ,m ^L 'which does 

not depend on -dy, and therefore, can be computed for each 
\m\, |n|, \m'\ < L with the complexity 0{A). Using Gm,n,m', 
Am,nAv) can be computed for all \m\, \n\ < L and for all 
sampling points dy using (21) in O(L^). Once Am,nAv) 
is computed, the signal f{(py,dy,ujyj), sampled over either 
of the proposed sampling scheme, can be obtained with the 
complexity O(L^logL) by employing FFT along m and n 
for each dy, as given in (10). The overall complexity of the 
resulting algorithm to compute IFT is therefore 0{A). 

F Computation of Wigner-d functions 

In the computation FT and IFT, we need to address the 
computation of the fixed-angle Wigner-d functions = 

d^ ^(7r/2), which are required to be computed for all £ < L, 
and all \m\, \n\ < 1. Let denotes a matrix of size (fll -f 
1) X {21 -f 1) with entries A^ „ for |m|, \n\ < i, which is 
computed for each f = 0, 1 using the recursion 

relation proposed in [25] that recursively computes A from 
jjjg computational complexity of 0{Lf). Since £ 
is of order L, the total complexity to compute Af„ „ for all 
Q < £ < L and |m|, \n\ < £ is 0{Lf). In the computation 
of both FT and IFT, this recursion is useful as A^ matrices 
are computed recursively for £ = 0, 1, ..., L — 1. In order 
to speed up the computation, we also employ the symmetry 
relations of Wigner-d functions [22], which facilitate us in a 
way that we are only required to compute (£ -f 1)^ entries 
of the matrix A^. The remaining of (2£ -f 1)^ entries can be 
inferred using the symmetry relations. 



Band-limit, L 

Fig. 1: Accuracy analysis: plots of the maximum error i£„iax 
and the mean error i£mean, respectively given in (22) and (23), 
for band-limits in the range 2 < L < 128. 


IV. Accuracy Analysis 


In order to analyse the numerical accuracy (and stability) 
of the propose sampling scheme and the associated Fourier 
transforms, we carry out a following experiment: 1) obtain a 
band-limited test signal, denoted by /t G 7 £l, by generating its 
Fourier coefficients {ft)m n for 0 < £ < L, |m|, \n\ < £ with 
real and imaginary parts uniformly distributed in the interval 
[—1, 1], 2) apply the IFT to obtain the band-limited test signal 
ft € TLl over the samples in the proposed scheme ©l), 3) 
using the samples of the test signal, the FT is then computed 
which yields the reconstructed Fourier coefficients, denoted by 
(/r)L.ri> 4) record the maximum reconstruction error i£„iax 
and the mean reconstruction error Fimeaiu given by 


L'max — max |(/t)m,n (/r)m,n|) 

1 


L-l 




= E E l(/.)G-(/.)L 




( 22 ) 

(23) 


The implementation the FT and IFT is carried out in double 
precision arithmetic in MATLAB. For a particular band-limit 
in the range 2 < L < 128, we repeat the experiment for 5 
different test signals and obtain the average values of i£„iax 
and Ejnean, which are plotted in Fig. 1, that demonstrates that 
the IFT of any band-limited signal followed by the FT yields 
the same band-limited signal with errors on the order of the 
numerical precision. Hence, the proposed scheme supports 
exact computation of Fourier transforms on the rotation group. 


V. Conclusions 

We have developed a Gauss-Legendre quadrature based 
sampling scheme on the rotation group for the discretisation 
of band-limited signal defined on the rotation group such 
that the FT can be exactly computed from its samples. The 
proposed GL sampling scheme is asymptotically as efficient 
as the most efficient scheme developed recently. For the fast 
implementation of the FT and IFT, we have also developed 
fast algorithms of complexity 0{A). We have analysed the 
computation time and demonstrated the accuracy of the FT for 
the proposed sampling scheme up to band-limit L = 128. 
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